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Abstract 

The following article has been submitted to/accepted by the American Journal of Physics. After 
it is published, it will be found at http://scitation.aip.org/ajp/ 

We provide a semiclassical optics derivation of Einstein's rate equations (ERE) for a two-level 
system illuminated by a broadband light field, setting a limit for their validity that depends on 
the light spectral properties (namely on the height and width of its spectrum). Starting from the 
optical Bloch equations for individual atoms, the ensemble averaged atomic inversion is shown to 
follow ERE under two concurrent hypotheses: (i) the decorrelation of the inversion at a given time 
from the field at later times, and (ii) a Markov approximation owing to the short correlation time 
of the light field. The latter is then relaxed, leading to effective Bloch equations for the ensemble 
average in which the atomic polarization decay rate is increased by an amount equal to the width of 
the light spectrum, what allows its adiabatic elimination for large enough spectral width. Finally 
the use of a phase-diffusion model of light allows us to check all the results and hypotheses through 
numerical simulations of the corresponding stochastic differential equations. 
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I. INTRODUCTION 



Lorentz's equatioiP and Einstein's optical rate equation^ are well known, very valuable, 
and widely used heuristic models for the description of light-matter interaction. Both of 
them occupy an important place in the early teaching of quantum optics topics too, each 
of the models applying to different situations and allowing the study of different physical 
problems (refractive index, laser action). Although both models were formulated before the 
advent of quantum mechanics, they can be justified a posteriori within the framework of 
quantum optics theory in the appropriate limits.'^ 

Generically speaking, formal derivations of heuristic models from first principles are im- 
portant not only for aesthetic reasons and completeness arguments, but also for the clarifica- 
tion of their applicability domains. In this respect the situation of Lorentz's and Einstein's 
models is quite different: derivations of the former from the optical Bloch equations are 
easily found in textbooks,'^ but this is not the case for Einstein's rate equations (ERE for 
short in the following), what results most surprising given the paramount importance of 
Einstein's model. Of course this does not mean that the connection between the optical 
Bloch equations and ERE has not been considered, as several quantum optics textbooks 
discuss ERE and provide derivations of Einstein's A and B coefficients for spontaneous and 
stimulated processes,^ see Refs. 4-6. More general treatments can also be found in some 
textbooks such as those of Refs. 4, 8. But these treatments do not focus on the derivation of 
ERE (rather on those of Einstein's A and B coefficients) and, from our viewpoint, do not 
put enough emphasis on didactic aspects. We try to close this gap with the present article. 

Einstein proposed his celebrated optical rate equations under the assumption of a strongly 
incoherent radiation, and this is the strict meaning of ERE. However, similar rate equations 
apply when the atomic line-width (alternatively, the decay rate of the induced electric dipole) 
is much larger than the decay rates of the atomic levels (see, e.g., Ref. 3), which occurs, in 
particular, in many laser systems and this is the reason why this approach is followed in 
most laser textbooks. In some quantum optics textbooks such rate equations are derived 
from optical Bloch equations through the adiabatic elimination of the medium polarization 
and considering a fully coherent radiation field — unlike Einstein. These derivations lead 
to correct rate equations but it is evident that they apply to situations (strong atomic 
incoherence) that are far from the ones where ERE do (strong radiation incoherence). That 
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difference manifests not in the form of the equations (both are rate equations involving 
the populations of the two atomic levels) but in the expressions of the coefficients appearing 
therein: Only in the limit of strong radiation incoherence the coefficient governing stimulated 
processes in the rate equations is Einstein's B coefficient. 

Our treatment below shows, among other things, that the effect of radiation incoherence 
manifests as an increase in the decay rate of the ensemble averaged atomic dipole. Hence 
when that effective damping rate is large as compared to the population inversion decay 
rate (see below a more rigorous statement) the adiabatic elimination of the electric dipole 
is justified, no matter whether the light incoherence or the atomic incoherence or both are 
responsible for that largeness. That procedure leads to optical rate equations in which the 
radiation-matter coupling constant depends both on the field spectral width and on the 
dipole decay rate, which however do not play a symmetric role. This explains the difference 
between an adiabatic elimination based on a large atomic incoherence and that based on 
a large light incoherence: Only in the latter Einstein's B coefficient is obtained. We also 
notice that in common textbook derivations of the B coefficient a weak field is assumed,'^ 
a limitation absent in our derivation below. 

The rest of this article is organized as follows. In Section II we present ERE briefly. In 
Section III we introduce the optical Bloch equations for a set of atoms and reduce them 
to an integro-differential equation for the ensemble averaged atomic inversion. Then, upon 
applying a decorrelation approximation between light and atoms, the population inversion 
dynamics gets directly connected to the light spectrum. This decorrelation approximation, 
which should hold for incoherent light fields, will accompany us along all our derivations. 
From that integro-differential equation, in Section IV we derive ERE by assuming addition- 
ally a Markov approximation, which is justified in the limit of very broad spectra, much 
broader than the atomic linewidth (strong light incoherence). Still adopting the Markov 
approximation, in Section V we generalize our previous analysis by considering that either 
the light spectrum is broad or the atomic line is (as compared to the inversion relaxation 
rate), or both, what allows us to discuss the combined role of light incoherence and atomic 
incoherence. In the first case the usual ERE are retrieved with the correct expression for 
Einstein's B coefficient while in the second case one recovers the rate equations used, e.g., in 
laser modeling. A different approach is used then in Section VI, where we remove the Markov 
approximation and transform the original integro-differential equation into effective Bloch 
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equations. Importantly in such equations the effective ensemble- averaged atomic coherence 
is shown to display a decay rate which is equal to the sum of its bare decay rate and the 
width of the light spectrum. Such effective Bloch equations are analyzed in different limits, 
in particular when the adiabatic elimination of the effective atomic coherence is in order. 
This leads to an alternative derivation of ERE, as well as to generalized rate equations, and 
allows setting the conditions under which rate equations (ERE in particular) actually hold. 
In Section VII we study numerically a particularly simple case (that of a light field having 
only phase noise), which enables us to show the validity of the statistical decorrelation as- 
sumption used in all the previous derivations and to discuss a number of questions. Finally 
the main conclusions of the work are given in Section VIII. 



II. EINSTEIN'S RATE EQUATIONS 

ERE describe the interaction of a broadband isotropic light field with a two-state atomic 
system. Einstein^ postulated three basic light-matter interaction processes (stimulated 
absorption and emission, and spontaneous emission), and established the rate equations 
governing the evolution of the populations of each of the atomic states. Denoting by Ni 
the population of the lower {i = 1) and upper {i = 2) atomic states, and assuming that 
iV = A^^i + A^2 is constant and large enough for individual absorptions and emissions only 
produce smooth temporal changes in Ni, the time evolution of the populations is given by 
(see Appendix I for a generalized form) 

dN2 



dt 

dNi dN2 



-AN2 + BW21 {Ni -N2), (la) 

(lb) 



dt dt ' 

where W21 is the spectral energy density of the light field at the atomic transition Bohr fre- 
quency U21, and A and B are Einstein's coefficients for spontaneous emission and stimulated 
processes, respectively, which neither depend on the field strength nor on time and verify 
A/B = hjjj2i/T^'^c' ■ In terms of the normalized population inversion n = {N2 — Ni) /N, 

hence — 1 < n < 1, Eqs. ([T]) have the simpler looking form 

dfi 

= -A{n + 1) -2BW2in. (2) 

Lit' 

We address the reader to Ref. 5 for a particularly didactic presentation and discussion of 
ERE. 
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III. BLOCH EQUATIONS FOR AN ENSEMBLE OF TWO LEVEL ATOMS 



Consider a generic light field whose electric component £ we write as 

£"'(r,t) = |[E(r,t)+E*(r,t)], (3) 

where E (r, t) is the field negative-frequency part (that containing terms oscillating as e"""^*, 
w > 0) interacting with a collection of identical two-level atoms or molecules (|2) and |1) 
will denote the atoms' excited and fundamental states) with Bohr frequency UJ21 and electric 
dipole matrix elements (2| /i |1) = (1| /i |2) = /iz, which have been taken to be real vectors, 
aligned parallel to the Cartesian z-axis, without loss of generality. Working in the Dirac 
picture, and after performing the rotating-wave approximation, the semiclassical optical 
Bloch equations for an individual atom, labeled by a and located at Fq,, can be written ad^ 

-^ = -A{n^ + l)-i {nia^ - n^al) , (4a) 

— = -7_l(Tq - -ilaUa, (4b) 

where Ua = P22 ~Pii (Xa = exp {iuj2it) denote, respectively, the population inversion 
and the slowly varying atomic coherence of atom a described by its density matrix p^"\ and 

^^a(t) = ^i5;,(r„,t)e''^-*, (5) 

is the complex Rabi frequency of the light field at the location of atom a, with E'^ = z • E. 

The effect of spontaneous emission has been phenomenologically included through the 
damping terms, as standard semiclassical theory cannot describe this process:^ We assume 
a decay rate A for the population inversion which implies that the decay rate of the atomic 
coherence should be \A as we are assuming that only the upper atomic state is affected by 
spontaneous emissionP However we shall attribute to a a a decay rate 7_l = \A + T'^'^ which 
includes an additional decay rate T^'^ describing the effects of dephasing collisions (those 
affecting the atomic coherence but not the population inversion). Here we do not consider 
radiative collisions (which affect both rta and o"q ) in order to keep the problem simpler, but 
they can be easily included (see Appendix I for a brief discussion on this). 

In order to cast our problem in a way similar to ERE ([2]), which involves just population 
inversions, we first eliminate the atomic coherence from (|4]) by integrating formally Eq. ( |4b[ ), 

a. (t) = -l f dt'n^ it') it') e-^^ (6) 
^ Jo 



where a transient term, cTq (0) exp (— 7j_t), has been dropped (alternatively one can take 
o"a (0) = without loss of generality assuming that the interaction is turned on at that 



instant). Plugging this into Eq. (4a) we get 

— = -A [n^ + 1) 

-Re f <lt'Vl^{t)ni{t')n^{t')e-'<^^'-''\ (7) 
Jo 

which is an integro-differential equation for the evolution of the population inversion of atom 
a. 



A. Ensemble averaging 

Equation ([t]) rules the population inversion dynamics of a single atom. As we are inter- 
ested in the average evolution of the whole system — the ensemble of atoms — , which is the 
quantity described by ERE, we introduce the ensemble averaged inversion 

n it) = {n^ (t)) , (8) 

where averages are defined as 

1 ^ 



and compute its evolution equation from Eq. ([T]) as 

-A{n + 1)- f dt'K{t,t')e-^^^'-''\ 
Jo 



dn 
dt 



:iOa) 



K (t, t') = Re {n^ (t) (f) {t')) . (10b) 



Equations ( 10 ) describe the average dynamics of the system in an exact way, in the sense 
that no approximation has been done on the original Bloch equations in order to arrive at 
them. Clearly it is necessary to evaluate the correlation function K {t, t') in order to perform 



the time integral in Eq. (10a) and then arrive at a connection of Bloch equations and ERE. 



B. The decorrelation approximation 



A direct comparison between Eqs. ^ and ( [lO| reveals a number of important differences. 
A main one is that Eq. ^ has no memory (the time derivative of n at time t just depends 



on the value of n at the same time t) while in Eq. (10a) memory effects are present. We 



will consider this point later because a previous issue is that in (10a) it is not only n (the 



ensemble averaged inversion) that rules its evolution but the individual atomic inversions, 



Ha-, through the compound correlation kernel K. A first necessary condition for ( 10 ) to merge 



with ERE is then that it should be possible to decorrelate K {t^t') in (10b) as 



K (t, t') ^ Re (fi, {t) it')) (n, (t')) = C {t, t') n (f) , 



(11) 



with 

C{t,t') = Re{Ql{t)Q^{t')), (12) 

the field autocorrelation function. Indeed this is a very reasonable approximation as K {t, t') 
contains correlations between 71^ {t') and Qa {t > t'), while {t) is assumed to have a ran- 



dom character. In the following we adopt this decorrelation approximation (11), which will 
be checked numerically in Section |VII 



Using (11), the Rabi frequency definition (|5|, and recalling that we are assuming that 
radiation is isotropic and unpolarized, which allows writing 



(E,(r,t)E:(r,t')) = |(E(r,t)-E* (r,t')) 



(13) 



Eq. (10a) becomes 



X Re / dt'n (t') (E (r, t) ■ E* (r, t')) e 



-(7x-i(^2i)(t-t') 



(14) 



This form is actually very interesting as it allows making direct contact with the light spec- 



trum via the Wiener-Khintchine theorem as shown in Appendix II. Making use of Eq. (48b) 
in that Appendix, i.e., 

2 r+oo 



we get 



(E (t) ■ E* (t')) = - / duW (u) e'"(*'-*), 

^0 J-oo 



^("-^' ^ ('-''' -^^'"''^ 



(15) 



1 r+oo 

I{t- t') = -Re / dooW (oo) e*('^2i-'^)(*-*')^ 



(16a) 
(16b) 
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where a factor vr has been included in the prefactor of the integral in (16a) for convenience, 
W (w) denotes the spectral energy density of the light field and I (t — t') is a frequency- 
shifted field autocorrelation function defined in terms of its spectrum. This is the equation 
we analyze throughout the rest of this paper. 

IV. A FIRST DERIVATION OF ERE 



As we already commented in the previous Section ERE have no memory while Eq. (16a) 
has. But let us put ourselves under the conditions considered by Einstein: If the light 
spectrum is broad (a highly incoherent light field) the field autocorrelation function I {t — t') 



defined in (16b) will be a very sharp function around t' = t. To be more precise, if we denote 
by A the width of W (u), then I {t — t') will be effectively zero but for \t — t'\ < tc, where 
the coherence time tc ~ A~^, as follows from the properties of the Fourier transform. Then, 
if 7_L ^ A as well one can substitute n (t') — )■ n (t) and e"'^-^^*"*') — )■ 1 under the integral 



in (10a). This constitutes a Markov approximation indeed as with it one assumes that the 
lack of correlation in the field (which implies a large enough amount of randomness in its 
evolution) provokes the complete loss of memory of the inversion n{t). Under this "light 



incoherence" dominated scenario Eqs. ( 16 ) become, after performing the time integration, 

_ = _A (n + 1) - j duW {uj) (5« {uj - uj2i) , (17) 

where 5*^*-' (x) = is a Dirac delta like function: It is a peaked function around x = 0, 
has a height equal to t/vr, a width equal to 27r/t, and verifies dx^^*-* (x) = 1. Hence, as 
soon as t ^ A^^ (which is a very short time), W {u) can be picked out of the integral as 
W {(jj = UJ21) = W21, leading to 

cif) 

= -A(n + 1) -2BW2in, (18a) 

at 

which coincides with ERE ^ and gives the correct result for B (see Ref. 5). 

V. LIGHT INCOHERENCE VS ATOMIC INCOHERENCE. FROM ERE TO 
LASER RATE EQUATIONS 

In our previous derivation of ERE we have assumed that the light spectrum was suffi- 
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ciently broad as to bring n outside the integral in (16a) — the Markov approximation — and 
to make the replacement e"'''^''*^*'-* — t- 1. Nevertheless one can still adopt the Markov ap- 
proximation without imposing 6""^-^ — )■ 1, and this is what we face in this Section. We 
are hence considering the possibility that, either because the light spectrum is broad or the 
atomic line is, or both, (i.e., max (7_l, A) ^ A), the function I {t — t') e"'''^ under the in- 



tegral in (16a) is strongly peaked around t' = t, allowing again the substitution n {t') — )■ n (t) 



under that integral. Hence, under the Markov approximation Eqs. (16) become, after per- 
forming the time integration. 



dn 
dt 



-A (n + 1) - 2CBW2in, 



;i9a) 



where B reads as in (18b), and we defined the dimensionless coefficient 

du W(uj) ' 



C 



7_L 1^21 



-C 



OO - UJ21 
1± 



where 



ax) 



1 1 



(19b) 



(19c) 



TT 1 + 

This expression for C neglects transient terms proportional to e"'*'-^* in the final result, which 



is a good approximation as soon as 'j±t > 1. The Lorentzian C y Z^^^^ j is centered at = uj2i 
and has a width (FWHM) equal to 2'j±, hence it represents the shape of the atomic line: 
Coefficient ( is given by the convolution of the (normalized) radiation spectrum with the 
atomic absorption spectrum and then provides a measure of the strength of the interaction. 

Remarkably the Markov approximation (in conjunction with the decorrelation hypothesis) 
leads naturally to a rate equation description of the system dynamics. Nevertheless this rate 



equation (19a) is not ERE (18) in general because of the presence of coefficient only if 
= 1 ERE are obtained. Hence in order to know whether ERE describe the population 
dynamics or not it suffices to study the behavior of C, as we do in the following subsections. 



A. The light incoherence hmit 

In the case when the spectrum is very broad as compared to 7_l (A 3> 7_l: the "light 



incoherence" limit considered in the previous Section) the Lorentzian C in (19b) acts as an 
effective Dirac delta by selecting, from W (w), just the portion around u = (X'21. Then, if the 
spectrum is a smooth function of u one can substitute W {u) by W21 and, after performing 



the remaining integral in we get ( = 1, corresponding to ERE (18), in agreement with 
our first derivation. 



B. The atomic incoherence limit 

In the opposite limit, i.e., when 7^ ^ A, the atomic linewidth is much broader than the 
light spectrum; hence one should expect radiation to behave as effectively coherent. In this 



the portion around the peak of W. By analogy to the previous case we call this the "atomic 
incoherence" limit. Denoting by Uq the frequency at the peak of W, coefficient ( takes the 
following expression 



Note that this "atomic incoherence" limit corresponds to the usual case treated in many 
laser textbooks, where the polarization decay rate 7_l is assumed to be large, which allows 
the adiabatic elimination of the atomic coherence (see next Section for an in-depth study 
of this technique) . 

C. A bridge between both limits 

To conclude this Section we consider a situation where one can study, in a continuous 
fashion, the combined role of the light incoherence and the atomic incoherence. For this a 
specific form of the spectrum must be chosen. We consider the usual Lorentzian form 



where A represents the width (FWHM) of the light spectrum. Note that we are using a 
resonant spectrum; detuned cases can be treated along similar lines. The final result reads 



case it is W [u) that acts as an effective Dirac delta in (19b) by selecting, from C 





(20) 




(21) 



c 



A 



(22) 



A + 27x' 
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Note that < C ^ 1- This expression contains, as special cases, the hght incoherence 



hmit, in which C = 1 and ERE (18) are recovered, and the atomic incoherence hmit where 
( = all this in agreement with our previous discussion. We see then that there is a 
continuous transition from the regime where ERE apply, dominated by light incoherence, 
and that typical of laser physics, dominated by atomic incoherence. Let us insist in that in 
both cases, and in intermediate cases as well, the system dynamics is of rate equations type. 



see (19a), but only in the light incoherence limit true ERE are obtained. 



VI. A DERIVATION OF ERE BASED ON THE ADIABATIC ELIMINATION OF 
EFFECTIVE BLOCK EQUATIONS. VALIDITY LIMITS OF ERE 



In the previous Sections we have been able to provide a semiclassical optics derivation of 
ERE, with the correct expression for the B coefficient, under the decorrelation and Markov 
approximations, whenever the latter is due to a strong light incoherence. In this case both 
assumptions are closely related as both are consequences of the randomness of the incoherent 
radiation field. 

In this Section we give an alternative derivation of ERE by relaxing the Markov approx- 
imation, which will give us relevant information about the role of the field spectrum (height 
and width) and of the dipole relaxation on the population dynamics. Let us then return to 



the general case represented by Eq. (16), in which only the decorrelation approximation has 



been done. In order to obtain some general, analytic result, a choice must be made for the 
form of the spectrum and we use again a Lorentzian one. For the sake of simplicity we shall 
also assume that its central frequency is resonant with the atomic transition, although this 
is inessential and the derivations that follow are easily generalized to the non-resonant case. 



Using then (21) the intensity (16b) becomes / (r) = ^APV2iexp (— |r|) and Eq. (16) can 
be written as 

dfz 

^ = ^A{n + l)-2CBW,,q, (23) 

where 



g(t)=7f / dt'n(t')e-^"(*-*'\ 



7f 



7± + 



(24) 
(25) 



the B coefficient reads as in (18b), and ( is given in (22) 
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Notice that if n varies slowly during a time interval then n (t') can be picked out from 



the integral in (24) at t' = t (this is the Markov approximation we made in the previous 
Section) and, after performing the remaining integral, q{t) = n (t) once 'yf^t ^ 1, leading 



to ERE ( 18 ) but with a modified coefficient for the stimulated processes when ( differs from 



unity as already discussed. 



The integro-differential Eq. (23) can be easily transformed into a pair of coupled linear 
differential equations by taking the time derivative of (24), 

df? 

— = -A{n + l)- 2CBW2iq, (26a) 
^ = 7f {n - q) , (26b) 

which are effective Block equations with q playing the role of a kind of normalized average 
medium polarization. Let us recall that this set of equations is exact (for a Lorentzian 



spectrum) but for the decorrelation assumption (11 ), which should hold under a wide variety 



of conditions. Furthermore notice that for A — 0, i.e., for a fully coherent field hence 



characterized by a constant Rabi frequency VLq (real without loss of generality), Eqs. (26) 
are equivalent to Eqs. ^ (for Via (t) = Va) upon identifying q with iAa/^lo and W21 
with Oq/A, as it must be.'^ 

We see that by assuming the decorrelation hypothesis and by taking a Lorentzian spec- 
trum for the radiation field, the average Bloch equations for the atom gas can be reduced 
to a pair of effective Bloch equations in which the effective medium polarization decays at 
a rate 7^^ = 7_l + A/2 = {A + A) /2 + r<^=, where A is the FWHM of the light spectrum. 
Within the range of validity of the assumptions the above amounts to say that the incoher- 
ence of the radiation field is in some way transferred to the average medium polarization, 
manifesting as an increase in its decay rate similar to the effect of non-radiative collisions, 
which are accounted for by F'^'^. Notice however that A and F'^'^ do not enter symmetrically 



in the equation of evolution of n, because of coefficient C, see (22). 



A. Adiabatic elimination of the effective coherence 



When 7^^ is large as compared with A — which is the situation we are interested in — (see 
below for a more rigorous bound) the effective coherence q can be adiabatically eliminated 



from the effective Bloch equations. We first integrate formally Eq. (26b) thus recovering 

12 



(24). Repeatedly integrating (24) by parts we easily get 



d2 



(27) 



and we see that for large enough we can approximate q{t) ^ n (t) (notice that this is the 



same as taking dq/dt = in Eqs. (26): the usual adiabatic elimination procedure). Using 



this result in the equation for n (t) we retrieve (19a), which was obtained under the Markov 
approximation in that Section. We see then that the latter and the adiabatic elimination 
of the effective coherence lead to the same result, as expected, thus completing our second 
derivation. Notice that, as already commented, we obtain a modified coefficient for the 
stimulated processes as in general C 7^ 1- We shall come back to this difference later. 

One of the virtues of this procedure is that it provides us with a simple tool to determine 
the conditions under which ERE (or rate equations in general) are correct, i.e., the conditions 
under which the above adiabatic elimination holds. 



We can estimate how large 7^^ must be by comparing the first two terms in (27). The 
necessary condition is |dn/dt| ^ 7^^ \n\, which using (26) leads to 



7f > A, BW21 < 



7f 
2C 



(28) 



for rate equations to be vahd. Using 7^^ = 7_|_ + A/2 = (A + A) /2 + r the above conditions 



read 



A + 2r^" > A 



BW21 < 



{A + A + 2T'^^y 
4A 



(29a) 
(29b) 



The second condition above sets an upper limit on the field energy density which is 
not usually stressed. The first condition implies that the adiabatic elimination is correct 
independently on which of the two quantities A or V^"^ is the larger one: The condition is 
that the effective coherence decay rate be large and not whether this is due to dephasing 
collisions or to light incoherence, in agreement with our initial discussions. 

However, importantly, the result is not the same for large A as for large T'^'^: In the 



limit of large light incoherence A ^ A, 2r'^^ one has C = 1 ^ind Eq. (19a) exactly coincides 
with ERE, whilst in the limit of large dephasing collisions rate, T'^'^ ^ A, A, one has ( ~ 
A/(2r^'=) < 1. This difference can be rephrased in the following way: For large A the 



13 



coefficient for stimulated processes is Einstein's B, whilst for large T the coefficient is not 
B but B'^^ = (B = A/ (2r*^'^) B (see endnote 10). This is the essential difference: Einstein's 
B coefficient appears only when the field is sufficiently incoherent and not when dephasing 
collisions are large. This is so because although dephasing collisions and light incoherence 
enter symmetrically in the equation for the effective atomic coherence q, through the effective 
damping rate 7, they don't in the equation for the inversion precisely because of the form 
of coefficient (. 

We note here that in most textbooks rate equations are derived from optical Bloch equa- 
tions ^ through the adiabatic elimination of the medium polarization without performing 
any ensemble averaging that accounts for the light incoherence. In this case a constant Rabi 
frequency is usually assumed, corresponding to coherent radiation, and the atomic dipole 
relaxation rate 7_l is assumed to be much larger than the population one because of the 
existence of frequent dephasing collisions. As we have seen this is a correct and legitimate 
way to derive rate equations, but it is not the right way for deriving ERE: The price paid 
with this simplified presentation is an incorrect expression for the B coefficient. 



B. Comparison of rate equations and effective Bloch equations solutions 



Compared to ERE (18) the system (26) has an extra equation that allows for a richer dy- 



namics. We shall now compare the solutions of both models in two different time regimes to 
obtain a better estimate of the necessary conditions for ERE be valid than that of inequality 

m. 



1. The short time limit 



For short times after the illumination has been switched on, at t = 0, the predictions of 



both models given by ERE (19a) and the effective Bloch Eqs. (26) differ 



Assuming n (0) = —1 and g (0) = we easily obtain 



nERE it) ^ -1 + 2(51^21^ + O {e) 
nsioch (t) ^ -1 + ]^CBW2it^. 



(30a) 
(30b) 



This means that ERE are overlooking the dynamics at the initial times, see Fig. 1(b), in a way 
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similar to what happens with the apphcation of Fermi's golden rule to the photo-ionization 
problem-fHl A way to cure the problem consists in taking n{t') out of the integral in (24) 
when considering the strongly incoherent limit (i.e., large 7^^), and retaining the exact value 
of the integral, i.e., approximate q (t) by ^1 — e^'^j"*^ n (t). With this approximation we get 
from (23) the modified ERE, 

dn 
dt 



-A (n + l)- 2CBW2in (t) (^1 - e"^!"*) , (31) 



that predicts a short time evolution as in (30b), see Fig. 1(b). We can now see clearly that 
the linear time dependence predicted by ERE ( 19a[ ) at short times is an artifact as we are 
using an approximate equation in a region {t < (7!^) ^) where it is not supposed to be 
valid. 



2. The long time limit 



At long times both the rate equations (19a) — ERE (18) is special case — and the effective 



Bloch Eqs. (26) reach the same steady state. 



n 00 



A 

~ A + 2CBW21' 

A stability analysis of ERE provides further insight in how this steady state is approached: 
It is performed by considering a situation in which n (t) is close to n (00) and seeing how 



the increment 5n (t) = n{t) — n (00) evolves according to (19a). One trivially gets 

d 



dt 



5n = XereSu, Aere = - (A + 2(^1^21) 



(33) 



which leads to a monotonous evolution 6n oc exp (Aere^) in which 6n decreases in time. 



The situation is slightly more involved in the effective Bloch equations case (26). Like 



before, we introduce the two increments 6n (t) = n(t) — n (00) and 6q (t) = q (t) — q (00) 
[Notice that q (00) = n (00)] and obtain 

-A -2CBW2i\ I 6n 




>,cff 



-if J \Sq 



(34) 



which leads to an evolution of the form exp {X±t) with 

- (7f + A)±Vn 



n 



-^-A)- 87f C5W^2i = ^ - 



(35a) 
(35b) 
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where in the last expression we used (22). In the strongly incoherent limit defined by 
7^*^ ^ A, 8BW21 the above eigenvalues take the simple form A+ — t- Aere, and A_ — )► —7^^: 
The large and negative eigenvalue A_ is responsible for a fast evolution that equalizes Sq 
and 6n, and from then on the system evolves as governed by ERE. 

For smaller incoherence, however, X± may become complex, thus signalling Rabi oscilla- 
tions with angular frequency equal to ^^/—7l, which requires TZ < 0, i.e., 

(A + 2r'i=-A)^ A , , 

BW21 > 7 —. 36) 

16A A»A,2rdc 16 ^ ^ 

Hence there are relaxation oscillations in the approach to steady state whenever the light 



spectral energy density is large enough, in the sense of Eq. (36), in stark contrast to the ERE 
prediction. Nonwithstanding, the observability of the oscillations must be examined because 



these are damped oscillations according to Eq. (35): In order that relaxation oscillations are 
present their frequency should be larger than, say, half their damping rate. It is easy to check 



that Eq. (36) already implies that condition, hence we conclude that relaxation oscillations 



will occur whenever Eq. (36) is fulfilled. 

We conclude that there will be no qualitative differences between the predictions of rate 
equations and effective Bloch equations, i.e., that rate equations are correct, whenever 

(A + 2rd=-A)' ^ ^ 

BW21 < ^ ^ ^, (37) 



which clarifies the meaning of symbol ^ in inequality (29). Of course, that these equations 
actually correspond to ERE also requires the light incoherence limit A ^ A, 2r'^^, in which 
case the condition becomes BW21 < ^■ 



VII. STOCHASTIC SIMULATION OF A FIELD WITH PHASE NOISE 

In the previous Sections we have derived ERE (and other rate equations) by making use 
of two hypotheses: the decorrelation between the field and the inversion and the Markov 
approximation. As the decorrelation hypothesis is crucial for all subsequent derivations 
and although reasonable because both Qa (t) and n (t) are random variables, we verify its 
validity for didactic purposes. For this aim we consider a light field having only phase 
noise and a tunable degree of incoherence, and take F'^'^ = in order to concentrate on the 
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role of light incoherence. We write the Rabi frequency at the location of atom a (|5]) as 
Qa {t) = fiocxp [—icpa (t)], with Qq = and the phase follows the evolution equation 

d0a 



dt 



(38) 



being A a diffusion coefficient that measures the degree of incoherence of the light field and 
(t) a random function without correlation between atoms. A simple case in which this is 
verified is when (^)}^=i are different realizations of white Gaussian noise ^ (t) which has 
a mean of (^^ (t)) = and correlatiori^ 



(39) 



Note that the averaging operator (|9j) has the same effect as the stochastic averaging because 
we are identifying each atom with a single realization of the problem, and for this reason we 
keep the same symbol, namely "()", for the averaging in both pictures. The problem can 
thus be described by a set of stochastic differential equations (SDEs), comprising the Bloch 
Eqs. and Eq. (g: 

dn 

dt 
da 

dt 

d0 

dt 



-A (n + l)- iQo [ae''*' - a*e-'^) 



(40a) 
(40b) 
(40c) 



Note that we have omitted the atomic subscript a as Eqs. (40) are interpreted as SDEs 



The field model used here represents laser radiation of finite linewidth.'^ Indeed the 
spectral energy density W (cj) of such field is given by, see Appendix [X| 



—Re 
B 



(/gi[</>(t+r)-(A(t)]\^ gi(a;2i-a;)r^ 



(41) 



where Eqs. (|5]), (13), (48a]) and (TSbj) have been used. As is a Wiener process, see Eq. (40c), 
its stochastic average is (exp [i [0 (t + r) — (t)])) = exp (— |A |r|) , see Ref. 9, and therewith 
we get 



(A/2)^ 
Afi (A/2)V(u;-a;2i) 



2 ' 



(42) 



which is a Lorentzian spectrum, as (21), centered at a; = uo2\ and having a width (FWHM) 
equal to A. In this case 

BW^x = f^o/A. (43) 
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Clearly A controls the degree of incoherence as anticipated. We note that ERE ^ in this 
case, according to (43), reads 



dn 
dt 



-A{n + 1) -2^n. 



(44) 



We integrated numerically Eqs. (40) by means of a fixed-step midpoint rule algorithm-fl^ 
The Gaussian noise ^ (t) was generated with the Box-MuUer method.'^ Depending on A, 
we averaged our results over a number of = 10^—10® realizations — or atoms — as many as 
necessary in order to obtain smooth results. We used different values of Qq and A and in 



all cases we fully confirmed the validity of assumption (11). We also compared the results 



of the full time evolution as given by the effective Bloch Eqs. (26) derived in the previous 



Section and Eqs. (40), as shown in Fig. 1 where the ensemble averaged population inversion 
n is represented as a function of the dimensionless time At for three values of the Rabi 
frequency. No difference can be found between the two predictions, which again confirms 



the validity of the decorrelation hypothesis ( 11 ). Hence we arrive at the conclusion that for a 



field whose incoherence is solely due to phase noise, Eq. (38), the effective Bloch model (26) 



is exact, and hence also that ERE (44) are exact for large enough A, roughly for A > 4f2o, A, 



see Eq. (37). 

We also paid attention to the influence of the number of stochastic trajectories A^. In 



Fig. 2 we represent the results of the numerical integration of Eqs. (40) for two different sets 
of parameters corresponding to cases with and without relaxation oscillations (see caption) 
using different values of (namely = 1, 10, 100, and 1000). The first feature to be noticed, 
see Fig. 2(a), is that for A^ = 1 and 10 the trajectories exhibit noisy Rabi oscillations which 
disappear for larger A^. In Fig. 2(b) the Rabi oscillations manifest as relaxation oscillations, 
but a closer look reveals that also in this case the individual Rabi oscillations manifest 
well beyond the disappearance of the relaxation oscillations. In other words, the individual 
atomic behavior exhibits Rabi oscillation, and it is the averaging that removes them. 

It is remarkable how the smooth averaged trajectories of Fig. 1 are approached as A^ grows, 
and it is interesting to note that for small A^ some oscillations are seen that disappear 
for larger values of A^. These plots suggest that A^ = 1000 is a large enough number of 



trajectories (or atoms) for the effective Bloch model (26) or ERE (44 ), to be valid, given 
the chosen parameters. This is more clearly seen in Fig. 3 where we represent the steady 
state reached by the system after a long enough transient, as well as its uncertainty, as a 
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function of N. Notice that for N > 10^ the uncertainty is almost neghgible. 

VIII. CONCLUSIONS 

We have provided a straight derivation of Einstein's rate equations (ERE) from the 
semiclassical optical Bloch equation for an ensemble of closed two-level atoms or molecules. 
The derivation has been done by assuming the statistical decorrelation between the inversion 
and field correlation function, and a Markov approximation owed to the (assumed) broad 
light spectrum. As well connections between ERE and usual laser rate equations have been 
considered, showing that leading to ERE which then have been analyzed in both the hmit 
where light incoherence dominates, resulting in ERE with the correct Einstein B coefficient, 
and the limit where atomic incoherence dominates, leading to laser rate equations. In 
the second derivation, the Markov approximation was replaced by the assumption of a 
Lorentzian spectrum in the radiation field. This second derivation led to a set of effective 
Bloch equations that contain information about the radiation spectrum whose bandwidth 
appears as an increase in the effective coherence decay rate. Then, for large enough spectral 
width, ERE are derived by adiabatically ehminating the effective coherence. Through the 
analysis and comparison of the solutions of both the ERE and effective Bloch models we have 
derived the conditions under which the former are applicable. We have discussed the subtle 
difference that exists between an adiabatic elimination based on large A (large spectral 
width) that led to ERE and provided the correct expression for Einstein's B coefficient, and 
an adiabatic elimination based on large 7_l (large atomic coherence decay rate) which leads 
to correct laser rate equations but does not provide a correct B coefficient. In other words: 
Einstein's B coefficient can only be correctly derived for large spectral width. We think 
that this is an important issue from the conceptual and pedagogical points of view. As well, 
upper bounds on the field strength have been derived which should be fulfilled in order that 
a rate equation description is valid. 

Finally, we have studied numerically the decorrelation hypothesis and checked the differ- 
ent predictions in the special case of a field having only phase noise. We think our derivations 
will help students in understanding more clearly how the ERE model can be justified and 
under which conditions it can be applied. 

We gratefully acknowledge help from Maria Gracia Ochoa with some numerical simula- 
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IX. APPENDIX I: INCLUSION OF RADIATIVE COLLISIONS 

Consider that besides spontaneous emission from the upper to the lower atomic state we 
also take into account the existence of atomic collisions able of forcing atomic transitions 
(usually referred to as radiative collisions). If we denote by 7jj the rate at which these 
collisions transfer population from level i to level j {'^ij depending on the temperature and 
on the collision cross-section of the atoms forming the gas), we can rewrite Eqs. ([T]) as 

-^ = -{A + 721) N2 + + BW21 (iVi - N2) , (45) 

and dA''i/dt = —(iN2/(it. Equation (|2| reads now 

dfi 

— = -7|j [n - neq) - 251^21^, (46a) 

7ll = A + 721 + 712, (46b) 

neq = -l + 2^, (46c) 

711 

where rigq is the equilibrium inversion in absence of radiation. 

By imposing that at thermal equilibrium W21 be given by Planck's formula^ it is easy 
to see that 712 = 721 exp (— /ia;2i/A;B^); with Boltzmann's constant and T the absolute 
temperature. This relation between 721 and 712 means that, in thermal equilibrium, for each 
collision induced atomic excitation there must be a corresponding deexcitation. 

We can also add radiative collisions in the optical Bloch equations in a similar way. Now 
we must writ^* 

dn 

= --yii (n„ - Ueq) - i (fi*cr„ - QaCr*J , (47a) 

^ = -l±0-a - ^^aria, (47b) 

with 7_|_ = 7||/2 + F'^"^. Notice that, in general, 712 ~ is a good approximation. 



Let us insist in that the above Bloch equations, and also Eq. (46), apply to a closed two 
level atomic system. If the system is assumed to be open (i.e., if relaxation processes can 
bring the atom into atomic states different from the two that interact with the light field) 
the relaxation terms need to be appropriately generalized (see, e.g., Ref. 3). 
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By following exactly the same lines we have followed in this article, the generalized 
Einstein Eq. (46) is derived from Eqs. (47). 



X. APPENDIX II: THE WIENER-KHINTCHINE THEOREM AND ITS APPLI- 
CATION TO EM RADIATION 

The Wiener-Khintchine theorem relates the spectral energy density with the field's au- 
tocorrelation function (see, e.g., Ref. 5). It can be put in either form 

p+oo 

W{io) = — dr(E(t) ■E*(t + r))e~*"", (48a) 

r) /•+00 

(E(t) ■E*(t + r)) = - / duW{u)e'''\ (48b) 

and this Appendix II is devoted to demonstrate this. The derivation bases on considering 
the e.m. field defined inside a fictitious volume (a cube of side L) with periodic boundary 
conditions allowing for the existence of traveling waves as in free space, and then letting 
L — )■ oo, which allows passing to the continuum. 

The electric and magnetic fields in such case can be written in their more general form 

as 

X(r,t) = J2Yl '^'n.ae*^'^""-''"*^ + ex., (49a) 

n (T 

where X = £,B, 

£n,a = L-^'^e, (n) (49b) 

-Bn.a = {cKY^ kn X £^^„, (49c) 

n G Z^, o" = 1, 2, kn = ^n, = |kn|, the polarization unit vectors (which we choose to be 
real -linear polarization basis- without loss of generality) verify e„ (n) • (n) = S^^y, and 
we take (— n) = e„ (n) by convention. First we compute the average e.m. energy density 
contained in the volume. 



u=--^{S'' + c' B ) = j-^ j d'r-^[£'' + c'B'']. (50) 



Upon using (49a) and taking into account that 

1 

[kn X (n)] • [k_n x e^/ (n)] = -klS^y, (51b) 



L3 
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one obtains, after little algebra, 

2£o5^5^L-='|A„,^|'. (52) 



u 



Finally passing to the continuum [d^/c = {2tx / Lf] and working in spherical coordinates 
[d'^/c = fc^d/cd^n] we get straightforwardly 

p+oo 

u= dooWiuj), (53a) 



JO 



J2 [ d^fi|A.(k)|L./,, (53b) 

„ Jin 



where du = cdk, and \Aa- (k)|^ = l^n.o-l^ such that n = k/ {2ti / L). 

Next we compute the correlation (E(r,t) •E*(r, t + r)), which following the previous 
lines can be written as 



1 2 iujT 



-1 p /• + 00 

(E (r, t) ■ W (r, t + r)) = — ^ V / d^fi / du;^;^ | A,, (k) 

zvr c ^ Jo 

Taking the Fourier transform of the above expression we get 

da;'(E(r,t).E*(r,t + r))e-'^=— 5^ / d'Q\AAk)\t^, 



which, compared with (53b), yields (48). □ 



XI. FIGURE CAPTIONS 
A. Figure 1 

Evolution of the ensemble averaged population inversion n as a function of the dimen- 
sionless time At. In (a) i^o = 4 and we used the three values of A indicated in the fig- 



ure. These results have been obtained by numerically integrating Eqs. (10) and the re- 



sults coincide exactly with those provided by the effective Bloch Eqs. (26). In (b) we used 



fio = V 11 {CBW21 = 2) and A = 5, and we have represented the predictions of the Bloch 



model (Eqs. (10), blue line), of ERE (Eqs. (44), red line), and of the modified ERE (Eqs. (31), 



brown line). The inset shows the different predictions for short times (see text). 
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B. Figure 2 



Evolution of the averaged inversion n obtained with the Bloch model Eqs. ( 10 ) for (a) 
Qq = 2 and A = 10, and (b) Qq = 6 and A = 1, for several values of the number of atoms 
A^. In (b) the trajectories for = 100 and = 10000 (not labeled for the sake of clarity) 
are so close each other that we plotted the former with dashed line in order to distinguish 
them. 



C. Figure 3 

Average population inversion ± its standard deviation after A^ trajectories for the steady 
state using (a) A = 10, fio = 2 (no transient oscillations) and (b) A = 1, fio = 6 (strong 
transient oscillations) . 
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